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A New Method for Recognizing Quadric Surfaces from Range Data and 
Its Application to Telerobotics and Automation 

(phase 1) 


by 
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Abstract 

In this phase of the proposed research a feature set of two-dimensional curves is 
obtained after intersecting symmetric objects like spheres, cones, cylinders, ellipsoids, 

paraboloids, and parallelepipeds with two planes. After determining the location and 

) 

orientation of the objects in space, these objects are aligned so as to lie on a plane 
parallel to a suitable coordinate system. These objects are then intersected with a hor- 
izontal and a vertical plane. Experiments were carried out with range images of 
sphere and cylinder. The 3-D discriminant approach has been utilized to recognize 
quadric surfaces made up of simulated data. Its application to real data has also been 
investigated. 
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INTRODUCTION 


Three-Dimensional object characterization and recognition is a primary task in the 
field of computer vision and automation. With the rapid improvement in quality and 
accuracy of several active sensors, digitized range data has proved to be the prime can- 
didate utilized for surface characterization. Several techniques have and are being 
developed, most of which rely heavily on differential geometry of smooth surfaces to 
describe 3-D objects. Although occlusion can be detected more easily in range images 
than in intensity based images, researchers have been trying hard to solve the problems 
of range image segmentation accurately and efficiently. Segmentation has been 
achieved for simulated and real range data for groups of symmetric objects of spheres, 
ellipsoids, cylinders and parallelepipeds. 

In this research we extend upon our earlier proposed research [1] wherein quadric 
surfaces are characterized by a feature set of two-dimensional curves. These curves are 
obtained after intersecting several 3-D surfaces with planes in various orientations. 
We put forward a new technique for determining the position and orientation of 
objects with respect to a fixed coordinate system. An efficient implementation of this 
algorithm will aid in obtaining a unique set of feature vectors for the various 3-D 
objects considered and thereby aid in the recognition process. The 3-D discriminant 
approach of classification and reduction of quadrics which was discussed in our earlier 
research [1,2], has been implemented and tested for simulated as well as synthetic 
range data. 

PREVIOUS WORK 

Since real range data obtained is not noise free, median filtering is performed [3] 
and the resulting images then undergo a series of curvature analysis tests to validate 
which image data can best be used for the recognition process. Best fit plots [2] also 
justifies the results obtained from curvature analysis as to which set of coefficients best 
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fits the range data. 

A quadric surface as shown in equation (1) is represented as a second degree 
polynomial, 

F(x,y,z) = ax 2 + by 2 + cz 2 + 2fyz + 2gzx + 2hxy + 2px + 2qy + 2rz + d = 0 (1) 

consisting of ten coefficients, which are obtained using the approach [4] developed by 
Groshong and Bilbro. 

Our proposed research has been applied earlier [1] to symmetric objects such as 
spheres, ellipsoids, cones, cylinders and parallelepiped. In the present research we 
extend it to various other quadric surfaces such as parabolic and hyperbolic cylinders, 
elliptic and parabolic paraboloids, and hyperboloids of one and two sheets. Each of 
these surfaces have been characterized separately in terms of its intersections with vari- 
ous planes. 

SURFACE CHARACTERIZATION 

We will consider the representation and intersections (by two planes) of the fol- 
lowing : parabolic and hyperbolic cylinders, elliptic and parabolic paraboloids, and 
hyperboloids of one and two sheets. Figure 1 illustrates each of the quadric surfaces 
to be used for the recognition process. 

Assume that all the objects are resting on planes parallel to the yz plane and that their 
axis of rotation is parallel to any one axis of the coordinate axis. Then in relation to 
equation (1) the product terms, namely yz, xz, and xy will be missing. We will refer 
to plane 1 as the plane that intersects the object parallel to the yz-axis, i.e., x = con- 
stant. Also let us refer to plane 2 as the plane that intersects the object parallel to the 
xz-axis, i.e., y = constant 

Let us initially consider the parabolic cylinder. This surface is symmetrical with 
respect to the xz plane the xy plane and the x-axis, and is a cylinder parallel to the y- 
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axis. Consider the cylinder lying on a plane parallel to the yz plane, thereby it is paral- 
lel to the y axis. The quadric representation without the product terms reduces to 

f(x,y,z) = az 2 + 2px + 2rz + d = 0 (2) 

Upon completing squares it reduces to 

(z + i) 2 + = 0 

a a 

only if d = - — . 

a 

Intersection of the parabolic cylinder with plane 1, i.e., x = k, 0 < k < 2p/ab, where b 
is any finite positive quantity, signifying the width of the base of the cylinder, yields 

(z + i ) 2 = 

a a 

which when solved results in a pair of parallel lines. 

Consider the intersection of the parabolic cylinder with the plane 2, i.e., y = k. 
Since equation (2) is independent of the variable y, the resultant curve intercepted is 
the same as equation (1), which is a equation of a parabola. 

The hyperbolic cylinder is considered next. This surface is symmetrical with 
respect to the coordinate planes and axes and any point on the z-axis. Since this 
cylinder is parallel to the x-axis, this variable will be missing in its representation. 
Equation (1) reduces to 


ay 2 + bz 2 + 2qy + 2rz + d = 0 (3) 

where ab < 0. Completing squares 


(y + -9-) 2 (z + -J-) 2 
a b 


a b 
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only if d = + — — 1. Also -1/b is a positive quantity. Intersection of the cylinder 

a b 

with plane 1, i.e., x = k, would generate a hyperbola. Since equation (3) is indepen- 
dent of the variable x, the curve intercepted is the one represented by equation (3). 
Intersection of the hyperbolic cylinder with plane 2, i.e., y = k would result in the 
equation 


(z+h 1 (k+^) 2 

b a 


which when solved results in a pair of straight lines. 

Let us consider the elliptic paraboloid next. This surface is symmetrical with 
respect to the xz plane, the xy plane and the x-axis. It is represented as 

ay 2 + bz 2 + 2px + 2qy + 2rz + d = 0 (4) 

Equation (4) upon completing squares, reduces to the form 


( y +^)2 (z+£) 2 

2 _ + !>_ + ^_=0 

1 1 1 


7 2 

only if d = — + 

b a 

316 sem * ma j° r and minor axes of the paraboloid, whereas l/2p is the 
height of the paraboloid. 

Consider the intersection of the elliptic paraboloid with the plane 1, i.e., x = k, 
where 0 < k < l/2r. Under these circumstances, the above equation results to 


(z + -£-) 2 (y + -^) 2 

-r^ + -r- 


6 


a 


(5) 



-k 


where — — is a positive quantity. Equation (5) is that of an ellipse. 

2p 

Consider the intersection of the surface with the plane 2, i.e., y = k, where — - \J — 

a V a 

< k < + a/X 

a V a 


Then 


(z » ^ 

b 

which is an equation of a parabola. 

The next 3-D surface considered is a hyperbolic paraboloid. Unlike the elliptic 
paraboloid, this object is symmetrical with respect to the xz plane, the xy plane and 
the x-axis. Its representation is as shown below : 

ay 2 + bz 2 + 2px + 2qy + 2rz + d = 0 
However in this case ab < 0. Upon completing squares we have 


(k + ^) 2 
a 


L 2p 


a 


(y + -^) 2 (z + ) 2 

a b 2p 


only if d = — + 


«*■<£ 


Intersecting the surface with the plane 1, i.e, x = k, results in 


( y + 1 ) 2 


(z + -J-) 2 
b 


zL 1 

a b 

where -1/a is a positive quantity. This equation is that of a hyperbola. In the case 


-k 

2p 


when x = k = 0, it will result in a pair of parallel lines which is a degenerate case of a 
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hyperbola. 

Consider the case when the object is intersected with plane 2, i.e., y = k, then 


(z % £) 2 

b 

b 


(k + - 1 ) 2 


a 


- 2px 


which is an equation of a parabola. 

The object considered next is a hyperboloid of one sheet. This surface is symmetrical 
with respect to each coordinate plane and axis and the origin. The equation of a 
hyperboloid lying in a plane parallel to yz plane is 


ax 2 + by 2 + bz 2 + 2px •+• 2qy + 2rz + d =0 

where we assumed the base to be a circular one and also ab < 0. 
Equation (6) upon completing squares results in 


( 6 ) 


-<x + -E-) 2 (y + h 1 (z+-j-> 2 
a + ±- + - b 


-1/a 


b 


b 


= 1 


F(x,y,z) = 


i_ j p 2 q 2 r 2 , 

where d = — f -- 1 — v 1 

abb 


If a < 0, i.e, b > 0, intersection of the hyperboloid with plane 1, i.e., x = k, where - p/a 


- *\J — < k < -p/a +*\j — , would result in 


( y +£) 2 (z + -^) 2 ( k + £ -) 2 

b + r*- = l + 


11 zL 

b b a 

where -1/a is a positive quantity. Equation (7) is that of a circle. 


(7) 


Intersection of the hyperboloid with plane 2, i.e., y = k, where - q/b - < k < 
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would generate 


(x + £) 2 (z +L ) 2 (k + -J-) 2 

a b b 1 

Zi 1 = zL 

abb 

where -1/a again is a positive quantity. This equation is that of a hyperbola. 

Unlike the hyperboloid of one sheet, the hyperboloid of two sheets consists of 
two separate pieces. Equation of a hyperboloid of two sheets lying on a plane parallel 
to the yz plane is 


ax 2 + by 2 + bz 2 + 2px + 2qy + 2rz + d =0 (8) 

where we assumed the base to be a circular one and also ab < 0. 

As in the case of the hyperboloid of one sheet upon completing squares, equation (8) 
results in 


-(X+-E-) (y + £) 2 ( z + t :) 2 

F(x,y,z) = 7 — — + : + t = -1 


-1/a 


b 


b 


2 2 j2 

where d = +-^-+— -+ 1. Intersection of the object with the plane 1, i.e, x 

abb 

where I k I > V-l/a, would result in 


= k. 


(y + £) 2 (z + £) 2 (k + £) 2 

b b , . a 


b b a 

where again -1/a is a positive quantity. This equation is that of a circle. 

However when I k I = V(- 1/a), the intersection will result in a point. 

Consider the case when the object is intersected with the plane 2, i.e., y = k, 
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where - q/b - < k < -q/b + 


Vf As in the case with the hyperboloid of one 


sheet, the intersection in this case results in 


(x + £) 2 (z + Lf (k + i)2 

-±~- f-^r +i 

abb 
which is an equation of a hyperbola. 

Table 1 summarizes the various curves (conics) evolved after intersecting each of 
the eleven surfaces with the two planes. 

TABLE 1 


OBJECT 

PLANE 1 : x = k 

PLANE 2 : y = k 

Ellipsoid 

ellipse 

circle 

Circular cylinder 

circle 

line 

Quadric cone 

circle 

hyperbola 

Sphere 

circle 

circle 

Parabolic cylinder 

line 

parabola 

Hyperbolic cylinder 

hyperbola 

line 

Elliptic paraboloid 

ellipse 

parabola 

Hyperbolic paraboloid 

hyperbola 

parabola,line 

Hyperboloid of one sheet 

circle 

hyperbola 

Hyperboloid of two sheets 

circle, point 

hyperbola 

Parallelepiped 

line 

line 




As seen from table 1 above, the quadric cone, the hyperboloid of one and two 
sheets all generate similar curves. Also the hyperbolic cylinder and the hyperbolic 
paraboloid generate the same curves after being intersected with planes 1 and 2. A 
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detailed analysis whereby the objects will be intersected with planes at various other 
orientation will solve the problem of objects sharing the same intersected curves. 

The surface characterization described in the previous few sections have been 
assuming that the objects are resting in a desired stable position, however in reality 
objects are disoriented in space. In the next few sections we put forward our new 
scheme which determines the orientation of the object in space and subsequently per- 
forms the desired rotations so as to place the objects in a stable position. 

DETERMINATION OF THE ROTATION MATRIX 

The determination of the location and orientation of a three-dimensional object is 
one of the central problems in computer vision applications. It is observed that most 
of the methods and techniques which look into this problem require considerable pre- 
processing such as detecting edges or junctions, fitting curves or surfaces to segmented 
images and computing high order features from the input images. Most of the 
methods are applicable to those objects which can be described by polyhedrons. The 
method we propose to investigate is based on analytical geometry whereby all the rota- 
tion parameters of any object placed in any orientation in space are evaluated sys- 
tematically. 

Let (x,y,z ) describe the coordinates of any point in our coordinate system. Con- 
sider a rotation of a about z axis, then the old coordinates in terms of the new one are 
represented as 


x = x'cosa + y'sina 

(9) 

y = -x'sina + y'cosa 

(10) 


i.e., the rotation matrix is 


Ra = 


cosa sina 0 
-sina cosa 0 
,0 0 1 


( 11 ) 




Next consider a rotation about the x' axis by an angle (3 of the same point. The resul- 
tant coordinates and the old coordinates are now related by the following equations 


z = z'cosP - y"sinP 

(12) 

y' = z'sinP + y"cosP 

(13) 


whereby the rotation matrix is 



1 

0 

0 


0 

cos(3 

-sin(3 


0 

sinji 

cosfr 


And finally consider a rotation about the y" axis by an angle y, then 


(14) 


z' = z"cosy - x"siny 


(15) 


x' = z"siny - x"cosy (16) 

The rotation matrix for the above process was 



cosy 0 siny 
0 1 0 


[-siny 0 cosy] 

Figure 2 illustrates the behavior of a particular point (x,y,z) with respect the various 
rotations described above. Using equations 11, 12, 13, 14, 15, and 16 solving for x, y, 
and z in terms of x", y", and z", yields 


x = x"(cosacosy - sinocsinPsiny) + y"sinacosP + z"(sinycosa + cosotsinctsinp) 
y = x"(-cosysina - sinasinPcosa) + y"cos(3cosa + z"(cosysinPcosa - sinysina) 
z = -x"sinycosP - y'sinp + z'cosycosP 

After substituting these new x, y, and z coordinates in the original quadric representa- 
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X 



x ” (d) 


Figures 2(a), (b), (c), and (d) refer to the coordinate system initially, 
after the first rotation, the second rotation , and subsequently the third 
rotation respectively. 
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tion, we get an entire set of new coefficients for x 2 , y 2 , z 2 , yz, xz, xy, x, y, and z. 

Our immediate aim is to eliminate the rotation terms xy, yz, and xz, namely F, G, 
and H which have been evaluated above. Let us initially list out the new coefficients 
for the product terms yz, xz and xy. 


2F' = j(bcos 2 a + asin 2 a + hsin2a - c)sin2p + (2gsina + 2fcosa)cos2pj cosy 
+ |((a - b)sin2a + 2hcos2a)cosP - (2gcosa - 2fsina)sinpj siny 

2G' = ja(cos 2 a - sin 2 asin 2 P) + b(sin 2 a - cos 2 asin 2 P) - ccos 2 p 

- (gsinasin2f3 + fcosasin2P) - hsin2a(l + sin 2 P) 

+ (a - b)sin2asinp + shsinpcos2a + cosP(2gcosa - 2fsina) cos2y 

2H' = [((a - b)sin2a + 2hcos2a)cosP - (2gcosa - 2fsina)sinpjcosy 

- j^(bcos 2 a + asin 2 a + hsin2a - c)sin2(3 + (2gsina + 2fcosa)cos2pj siny 


(17) 


(18) 


(19) 


Since our aim is to eliminate the rotation terms xy, yz, and xz, let’s now 
exclusively consider the coefficients of these rotation terms, namely F, G, and H which 
were evaluated above. In an iterative procedure we will be able to eliminate all of the 
product terms. For e.g., let’s consider we wish to eliminate the term xy. The, by a 
specific rotation of a about the z axis, we will be able to accomplish our goal. How- 
ever while carrying out this process, the orientation of the object about the two planes 
yz and zx, i.e., the angles the object made with these two planes have been disturbed. 
Now if we wish to eliminate the yz term, the object has to be rotated about the x axis 
by an angle beta. However in this instance, while carrying out the process, the already 
missing xy term reappears again, although the magnitude of its present orientation has 
been reduced. Hence by carrying out the above process in an iterative fashion, there 
comes an instance when all the coefficients of the product terms converge to zero. 

Let us consider the equations (17), (18), and (19) respectively. Let us eliminate 
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the coefficient h, i.e, the xy term in step 1. This can be accomplished by rotating the 
object about the z axis by an angle a, whereas (5=7=0. Under these circumstances the 
new coefficients look like as shown below 

2f n = 2gsina 1 + 2fcosa 1 (20) 

2g u = 2gcosa! - 2fsina l (21) 

2h n = (a - b)sin2aj + 2hcos2a 1 = 0 when cot2a! = b - — (22) 

2h 

As seen above the coefficient "h" has been forced to 0. The most significant bit of the 
subscript refers to the iteration number , whereas the least significant bit of the sub- 
script reflects the number of times the object has been rotated by a specific angle. In 
the above case the LSB of 1 refers to the first instance where the object has been 
rotated by an angle a. The remaining coefficients a, b, c, p, q, and r also reflect 
changes brought about by the above rotation. 

The new quadric equation now has a look as shown below: 

F(x,y,z) = a u x 2 + b^y 2 + c u z 2 + 2f n yz + 2p u x + 2q u y + 2r u z + d = 0 (23) 

Consider the second step wherein the coefficient corresponding to the yz term is 
forced to zero. In this particular case, the object has to be rotated by an angle (5 about 
the x axis, whereas a=y=0. Under these circumstances, the new rotation coefficients 
(signifying the product terms) becomes 


2f 12 = (b 12 - c 12 )sin2P 1 + 2^^082^! = 0 if 


0012 ^ = 



(24) 


2gi2 = Zgncos^! 


(25) 



( 26 ) 


2hj2 — — 2g^sinpj 

The new quadric equation as before looks like as shown below: 

F(x,y,z) = a 12 x 2 + b 12 y 2 + c 12 z 2 + 2f 12 yz + 2p 12 x + 2q 12 y + 2r 12 z + d = 0 (27) 

In the final step of the initial iteration, the coefficient corresponding to the xz term is 
forced to zero. In this case, the object is to be rotated by an angle y about the y-axis, 
whereas a=p=0. Under these circumstances, the new rotation coefficients become 

2f 13 = 2h 12 sinYi = -^gnsinf^sinYi (28) 

2gj3 =(a 13 - c 13) sin2Yi + (2g 11 cosa 1 -2f n sina 1) cosp 1 cos2Yi = 0 (29) 

Cj 2 — a I2 

if cot2Yi = g w 

2h 13 = 2h 12 cosYl = ^gusin^cosYl (30) 

Let’s now carefully analyze the coefficients of xy, yz, and zx obtained in the final step 
of the first iteration. Consider for instance the coefficient corresponding to the yz 
term. It is observed that while proceeding from one step to the other, the new 
coefficients are getting multiplied by the sine or cosine of the concerned angle. This 
implies that in every succeeding steps, these coefficients are decreasing in their magni- 
tude. To justify the above statement, let us now consider all the coefficients obtained 
in the second iteration. 

At the end of stage 1 of the second iteration, the rotation coefficients become 


2f 21 = 2f 13 cos02 = -2g u sin(3isinYiC0SO2 

(31) 

2g 2 i = -2f] 3 sina2 = 2g n sin|$lsinYisina2 

(32) 

bi 3 — a n 

2h 21 = 0 if cot2ct2 = h 13 

(33) 
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At the end of the second stage of the second iteration, the rotation coefficients become 


_ c 2 i - b^i 

2f a = 0 if C0t2p2 = - f 2 i 

(34) 

2g22 = 2g 11 sinP 1 sinYiSina 2 cos(3 2 

(35) 

2h 22 = -2g 11 sinP 1 sinY 1 sina 2 sin{3 2 

(36) 


At the end of the final stage of the second iteration, the rotation coefficients reduce to 


2fz3 = -2g 11 sinp 1 sinYiSinot 2 sinP 2 sinY 2 

(37) 

bii ~ 3i3 

2g23 = 0 if 0^20^ = h 13 

Jm r 

(38) 

2h 23 = ^gjjsinPiSinYisinc^sinJ^cosYz 

(39) 


where o^, P 2 , and y 2 are the respective rotation angles along the z, x and y axes in the 
second iteration. Hence it is observed that with each iteration the rotation coefficients 
get smaller and smaller in magnitude and eventually drop out. 

We are now in a position to formulate a rotation matrix whose elements 
correspond to the direction cosines of the x, y, and z axes of the rotated object. 


Rotation Matrix =R y RpR a 


(40) 


where 


Ra = 


cosa 

-sina 

0 


sina 0 
cosa 0 
0 1 



1 0 
0 cos(3 
0 -sinp 


0 

sin (3 
cos(5 


(41) 


(42) 
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(43) 


R,= 


cosy 0 siny 
0 1 0 
-siny 0 cosy. 


Subsequendy, 




cosacosy + sinasin^siny 
-cos^sina 

-sinycosa + cosysinasinJ3 


cosysina - sinysm^cosa 
cos^cosa 

-sinasiny-cosysin Pcosa 


where 


sinycosP 

cosP 

cosffcosy 


n n n 

a = 5X P = XPi» ^ Y = XYi- n corresponds to the iteration where all the rotation 

i=l i=l i=l 

terms go to zero. 

Once the rotation terms, i.e., xy, yz, and xz are eliminated, the 3-D surface has 
the representation of 


F(x,y,z) = Ax 2 + By 2 + Cz 2 + 2Px + 2Qy + 2Rz + D = 0 (44) 

where A, B, C, P, Q, and R are the coefficients evolved after the elimination of the 
rotation terms. As seen from above each of the product terms in the succeeding terms 
are getting multiplied by the sine or cosine of a concerned angle. This implies that in 
every succeeding steps, these coefficients are decreasing in their magnitude. Hence it 
remains to be seen that at which point, i.e, after how many iterations all the rotation 
coefficients converge to zero. 

We are aware that the intersection of a solid with a 2-D plane generates a curve. 
What type of curve it intersects, depends on the surface and the orientation of the 
plane of intersection. Since we have no idea about the type of surface, the only thing 
we can do is to experiment intersecting the surface with a series of planes. Does that 
mean we have to intersect the object with planes oriented at 1° apart from each other? 
This will mean we have to intersect the object 360 different times. 
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Our most important assignment will be to derive a consistent method for deter- 
mining the minimum number of planes that are necessary to intersect a given 3-D sur- 
face so that the generated conics uniquely characterizes the surface. This includes the 
derivation and formulation of the angular bounds for which a particular plane intersect- 
ing a surface results into the same 2-D curve. 

Consider the equation of an ellipsoid (spheroid) resting on a plane parallel to yz 
plane and its axis of revolution parallel to y-axis : 

The general representation reduces to the form 

F(x,y,z) = ax 2 + by 2 + az 2 + 2px + 2qy + 2rz + d = 0. (46) 

which further reduces to 


“ 

2 

* 

2 

- 

X + -2- 


y + S. 


Z + — 

a 


b 


a 




— + - 



a 


b 


a 


(47) 


2 2 j 2 

only if d = ^— + -^- + 1 and also a > 0, b > 0 . 

aba 

It should be noted that the coefficients a, b, p, q, r, and d are all known, and 

V? are the semi-major and minor axes of the ellipsoid, and [-p/a, -q/b, -r/a] are the 
coordinates of the center of the ellipsoid. 

Consider the intersection of the ellipsoid with plane 1, i.e., x = k, where 

=SL. JI <k< 3L + JI, 

a "a a “a 

then, 


(y + -J -) 2 (z + -) 2 

b a 

1 (ak + p ) 2 i 

1 

ak + p 

T 

i 

1 * 

■8 

JO 

a 



L 19 J 



( 48 ) 



which is the equation of an ellipse. 


Let’s now consider the intersection of the ellipsoid with plane 2, i.e., y=k, where 

then, 

(x + -2-) 2 (z + — ) 2 

T + T -1 = 0 (49) 

J_ _ (bk 4- q) 2 J_ _ (bk 4- q) 2 

a ab a ab 

which is the equation of a circle. 

To what extent the above mentioned two planes are able to distinguish the various sur- 
faces is to be investigated. 

Since different quadrics may generate similar curves, a detailed analysis whereby 
the objects will be intersected with planes at various other orientation will solve the 
problem of objects sharing the same intersected curves. Regions wherein a particular 
3-D surface yields a particular type of curve is to be observed and finally we hope to 
represent all the 3-D surfaces mentioned earlier by a unique five tuple vector. 

3-D DISCRIMINANT 


In the following few sections we will investigating in detail a 3-D approach of 
classification and reduction of quadrics [5], which looks into the various invariants of 
the quadratic form under translation and rotation of 3-D objects. 

The general equation of the second degree in three variables x, y, and z can be 
written in the form 

F(x,y,z) = ax 2 + by 2 + cz 2 + 2fyz + 2gzx + 2hxy + 2px + 2qy + 2rz + d = 0 
Associated with f(x,y,z) are two matrices "e" and "E", such as. 



e = 


a h g 
h b f 
J f c 


a h g p 
h b f q 
g f c r 

■P q r d. 


Let the determinant of "E" be denoted by A, and the determinant of "e" be 
denoted by "D". Also let the cofactors of each element of A be denoted by the 
corresponding capital letters. Based upon the fact that "E" is singular or non-singular, 
the 3-D surfaces can be singular or non-singular. Example of singular surfaces are 
ellipsoids, hyperboloids, paraboloids. The other quadrics are singular. 

It has been shown in [5] that I, J, D, and A, where I = a + b + c, J = a+p + y, 
D, and A the determinants, are invariant for any general coordinate transformation. 

J I 2 

Combining these four invariants a further set of three absolute invariants, I — , — , 

and is also obtained. 

A 

Based upon the above given set of invariants the following classification as sum- 
marized in table 2 is obtained. 
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TABLE 2 


Number 

Surface 

Rank 

Sign of A 

K’s same sign 

P3 

P4 

1 

Real ellipsoid 

3 

4 


2 

Hyperboloid of one sheet 

3 

4 

+ 

3 

Hyperboloid of two sheets 

3 

4 

” 

4 

Real quadric cone 

3 

3 


5 

Elliptic paraboloid 

2 

4 

“ 

6 

Hyperbolic paraboloid 

2 

4 

+ 

7 

Real elliptic cylinder 

2 

3 


8 

Hyperbolic cylinder 

2 

3 


9 

Parabolic cylinder 

1 

3 



Results obtained using the above procedure for simulated and real objects will be dis- 
cussed in the next section. 

EXPERIMENTAL RESULTS AND CONCLUSIONS 

For experimental purposes, range images of a sphere and a cylinder are con- 
sidered. Four types of data, namely, raw image of sphere, averaged image of sphere, 
raw image of cylinder and averaged image of cylinder were available to us. The raw 
data being noisy, had to be cleared using median filters. Mask sizes of 3 x 3 and 5 x 
5 were used to remove the noises and the effect of these filters on the range images 
were studied. Figure 3 shows the raw range image of the sphere. Figure 4 shows the 
raw image of a cylinder. Figures 5 and 6 are the segmented range images of the 
sphere and cylinder respectively. Effect of the 3x3 and 5x5 median filters on the 
segmented range images of the sphere and the cylinder can be seen in figures 7, 8, 9, 
and 10 respectively. 
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Figure 4. Raw range image of the cylinder. 
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Figure 5. Segmented raw range image of the sphere. 



Figure 6. Segmented raw range image of the cylinder. 
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Figure 1.3X3 median filtered image of the raw sphere. 



Figure 8. 5 X 5 median filtered image of the raw sphere. 
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Figure 9. 3 X 3 median filtered image of the raw cylinder. 



Figure 10. 5X5 median filtered image of the raw cylinder. 
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Curvature sign maps to check the integrity of the original as well as the processed 
images are pursued next. Not much changes were observed in the processed images. 
Figures 11, 12, 13, and 14 are the first and second derivatives w.r.t x and y axis for 
the image in figure 3. 

In all of these figures, the sign "+" is assigned to a particular pixel position if the 
magnitude of the derivative (first or second) of that pixel is greater than the magnitude 
of the derivative (first or second) of the pixel to its right. Similarly the sign is 
assigned to a particular pixel position if the magnitude of the derivative (first or 
second) of that pixel is less than the magnitude of the derivative (first or second) of 
the pixel to its right. In the case when the magnitudes of the derivatives (first or 
second) of either pixels is the same, the sign " " (blank) is assigned. 

Two sets of coefficients were generated for each set of the range data. While obtain- 
ing the best fit plot for the raw sphere, numeral 1 refers to the situation when the first 
set of coefficients of the raw image fits best, numeral 2 refers to the case when the 
second set of coefficients of the raw image fits best, numeral 3 refers to the case when 
the first set of coefficients of the 3 x 3 filtered image fits best, numeral 4 occurs when 
the second set of coefficients of the 3 x 3 median filtered image fits best, numeral 5 
occurs when the first set of coefficients for the 5 x 5 median filtered images best fits 
the data, and finally numeral 6 occurs in the plot when the second set of coefficients of 
the 5 x 5 median filtered image fits best to the data. 

Tables 3 and 4 tabulate the coefficients obtained for the various images of the raw 
range sphere. Figure 15 shows the best fit plot for these sets of coefficients. A similar 
procedure is carried over for the averaged range data of the sphere. 
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Figure 11. First derivative with respect to x-axis for the sphere raw 
image. 
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Figure 12. 
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First derivative with respect to the y-axis for the sphere raw image. 
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Figure 13. Second derivative with respect to the x-axis for the sphere raw image. 
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Figure 14. Second derivative with respect to the y-axis for the sphere raw image. 
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Figure 15. Best fit plot for the sphere raw image. 
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TABLE 3 


COMPARISON OF COEFFIC 

IENTS EVALUATED FOR 

THE ORIGINAL AND THE PROCE 

SSED IMAGES, SET 1 





COEFFICIENT 

RAW IMAGE 

3x3 FILTERED IMAGE 

3x3 FILTERED IMAGE 

2 

A. COBH 7 . OF X 

a01393 

ai»7 

0.1409 

B. COEFF. OF y 1 

-0.10112 

02433 

-01392 

C. COEFF. OF 

asm 

a7712 

06497 

F, COEFF. OF YZ 

08841 

•05509 

-03292 

0. COOT. OF XZ 

-03664 

04908 

•09748 

H, COEH 7 . OF XY 

-06073 

02561 

ai339 

P. COEFF. OF X 

aii36 

01857 

0.7426 

Q. COEFF. OF Y 

-0A968 

0.4079 

0.1361 

R, COEFF . OF Z 

-0J254 

43939# 

-08208 

D. CONSTANT 

aim 

03130 

02293 
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TABLE 4 


COMPARISON OF COEFFIC 
. — 

IENTS EVALUATED FOR 

THE ORIGINAL AND THE PROCE 

SSED IMAGES, SET 2 





COEFFICIENT 

RAW IMAGE 

3x3 FILTERED IMAGE 

3x3 FILTERED IMAGE 

2 

A. COEfF. OPX 

0.3108 

02027 

03030 


0.1298 

OI475 

04)392 

C, COEFP. OP 

0.6218 

0.7869 

06526 

F. COEFF. OF YZ 

-0.4998 

-03950 

•0.4487 

0. COEFP. OF XZ 

-0.8CQ2 

-06658 

•08376 

H, COEFP. OF XY 

0.3236 

01909 

02416 

P, COEFF. OF X 

0.4048 

033*4 

OA047 

Q, OOEff. OF Y 

0.2883 

02665 

02140 

R. COEFF. OF Z 

-0.6945 

-0.9176 

-0.7089 

D. CONSTANT 

0.1990 

02768 

0.1913 


In the case of the cylinders, as with the sphere, curvature analysis, the best fit 
analysis etc., are carried over to see which set of coefficients can be used for the 
recognition process. Coefficients generated for the various images of the cylinders can 
be seen in appendix A. 

The best fit plot for the raw sphere image prompts us to use the 3 X 3 filtered 
image of the second set of coefficients for the recognition process. Shown in tables 5 
and 6 are the various curves intercepted with the two planes for each of these images. 
Each of the objects, the raw and the processed images of the sphere and the averaged 
and the processed images of this sphere are intersected with two planes (one parallel to 
the yz-plane, the other parallel to the xz-plane). A decision on the curve being an 
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ellipse or a circle was made based upon the parity and the disparity of the x 2 , y 2 , and 
the z 2 coefficients. 

TABLE 5 


Sphere Images 

Raw 

3 x 3 filt. 

5x5 filt. 

plane 1, x = k 

Ellipse 

Ellipse 

Circle 

plane 2, y = k 

Ellipse 

Ellipse 

Circle 


TABLE 6 


Sphere Images 

Raw 

3 x 3 filt. 

5x5 filt. 

plane 1, x = k 

Ellipse 

Ellipse 

Circle 

plane 2, y = k 

Ellipse 

Ellipse 

Circle 


Similar experiments were carried over for the raw and the averaged images of the 
cylinder. In the case of the cylinder, since its orientation is unknown, we applied the 
rotation algorithm. This algorithm eliminated the product terms and then we con- 
ducted the same set of experiments as we had done before. Results for this procedure 
are listed in the tables 7 and 8. 
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TABLE 7 


Cylinder Images 

Raw 

3 x 3 tilt. 

5 x 5 tilt. 

plane 1, x = k 

Ellipse 

Line or Ellipse 

Line 

plane 2, y = k 

Ellipse 

Ellipse 

Ellipse 


TABLE 8 


Cylinder Images 

Raw 

3x3 tilt. 

5 x 5 tilt. 

plane 1, x = k 

Ellipse 

Line or Ellipse 

Line 

plane 2, y = k 

Ellipse 

Ellipse 

Ellipse 


As mentioned in the earlier sections, we had implemented the 3-D discriminant 
approach for several simulated as well as synthetic data. Results for the simulated 
data are listed in appendix B, and seem to be very effective as predicted by the theory. 
However, very unsatisfactory results were obtained while experimenting with real data. 
Appendix C contains the program listing for the various algorithms implemented. 
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The input file was "strawmedl . cod 
The output file is n strawmedl . coe 


The 

coef f 

of 

x-squared 

is 

-1.1279882E-02 

The 

coef f 

of 

y-squared 

is 

9 . 93 62820E-02 

The 

coef f 

of 

z-squared 

is 

0.9533964 

The 

coef f 

of 

yz 

is 

6 . 1207127E-02 

The 

coef f 

of 

zx 

is 

0.3908886 

The 

coef f 

of 

xy 

is 

7 . 4366599E-02 

The 

coef f 

of 

X 

is 

-0.4194884 

The 

coef f 

of 

y 

is 

-8 . 1451185E-02 

The 

coef f 

of 

z 

is 

-1.664232 

The 

constant 

d 

is 

0.750917.0 
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The input file was °strawmed2 . cod " 

The output file is ■ strawxned2 . coe " 

The coeff of x-squared is -6 . 1688922E-02 
The coeff of y-squared is 0.3596755 

The coeff of z-squared is 0.3823600 

The coeff of yz is -0.7434295 

The coeff of zx is -0.6636371 

The coeff of xy is 0.6694452 

The coeff of x is 0.8502086 

The coeff of y is 0.6296228 

The coeff of z is -0.6596514 

The constant d is 0.2292093 
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The input file was * STRMED3 1 . COD 
The output file is ■ STRMED3 1 . COE 


The 

coef f 

of 

x-squared 

is 

-0.5819576 

The 

coef f 

of 

y-sguared 

is 

-2 . 5062649E-02 

The 

coef f 

of 

z-squared 

is 

-0.4078553 

The 

coef f 

of 

yz 

is 

-9 . 12 91323E-02 

The 

coef f 

of 

zx 

is 

0.9860787 

The 

coef f 

of 

xy 

is 

8 . 9543425E-02 

The 

coef f 

of 

X 

is 

-0.3951663 

The 

coef f 

of 

y 

is 

4 . 5018956E-02 

The 

coef f 

of 

z 

is 

0.3026683 

The 

constant 

d 

is 

-4 . 4815123E-02 
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The input file was ■ STRMED3 2 . COD 
The output file is "STRMED32 .COE 
The coeff of x-squared is 0.1394030 

The coeff of y-squared is 0.2895423 

The coeff of z-squared is 0.4720904 

The coeff of yz is -0.6865588 

The coeff of zx is -0.7503323 

The coeff of xy is 0.5597892 

The coeff of x is 0.6822208 

The coeff of y is 0.5422743 

The coeff of z is -0.6979018 

The constant d is 0.2410266 
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The input file was "STRMED51 .COD 
The output file is ■ STRMED5 1 . COE 


The 

coef f 

of 

x-squared 

is 

6 . 2892713E-02 

The 

coef f 

of 

y-squared 

is 

5 . 5489447E-02 

The 

coef f 

of 

z-squared 

is 

0.9393368 

The 

coef f 

of 

yz 

is 

3 . 8853236E-02 

The 

coef f 

of 

zx 

is 

-0.4683146 

The 

coef f 

of 

xy 

is 

1 . 9889891E-02 

The 

coef f 

of 

X 

is 

0.3011477 

The 

coef f 

of 

y 

is 

-4 . 1780550E-02 

The 

coef f 

of 

z 

is 

-1.364659 

The 

constant 

d 

is 

0.5040251 
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The input file was * STRMED52 .COD ■ 

The output file is *STRMED52 .COE " 

The coeff of x-sguared is 5 .7214227E-02 

The coeff of y-sguared is 0.5990739 

The coeff of z-sguared is 0.4416272 

The coeff of yz is -0.8076080 

The coeff of zx is 0.4595302 

The coeff of xy is -0.1490175 

The coeff of x is -0.5915767 

The coeff of y is 1.089848 

The coeff of z is -1.019888 

The constant d is 0.6643057 
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The input file was ■ STAMED3 1 . COD 
The output file is ■ STAMED3 1 . COE 


The coeff of x-squared is 0.2759137 

The coeff of y-squared is 2 . 7527343E-02 


The 

coeff 

of 

z-squared 

is 

0.7029013 

The 

coeff 

of 

yz 

is 

0.1449835 

The 

coeff 

of 

zx 

is 

-0.9098228 

The 

coeff 

of 

xy 

is 

-9 . 6383080E-02 

The 

coeff 

of 

X 

is 

0.5634921 

The 

coeff 

of 

y 

is 

-8.9731783E-02 

The 

coeff 

of 

z 

is 

-0.8506840 

The 

constant 

d 

is 

0.2536311 
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The input file was "STAMED32 .COD 
The output file is " STAMED3 2 . COE 


The 

coef f 

of 

x- squared 

is 

0.6636790 

The 

coef f 

of 

y-squared 

is 

2.0903163E-02 

The 

coef f 

of 

z-squared 

is 

-9 . 2439458E-02 

The 

coef f 

of 

yz 

is 

-2 . 1956543E-02 

The 

coef f 

of 

zx 

is 

-0.7604508 

The 

coef f 

of 

xy 

is 

0.7227234 

The 

coef f 

of 

X 

is 

0.4242096 

The 

coef f 

of 

y 

is 

-0.2155603 

The 

coef f 

of 

z 

is 

0.3749632 

The 

constant 

d 

is 

-0.2534057 
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The input file was "STAMED51 .COD 
The output file is "STAMED51 .COE 


The 

coef f 

of 

x-squared 

is 

0.1115851 

The 

coef f 

of 

y-squared 

is 

3 . 13 68352E-02 

The 

coef f 

of 

z-squared 

is 

0.8936580 

The 

coef f 

of 

yz 

is 

0.1347357 

The 

coef f 

of 

zx 

is 

-0.5961419 

The 

coef f 

of 

xy 

is 

-4 . 8396215E-02 

The 

coef f 

of 

X 

is 

0.4117958 

The 

coef f 

of 

y 

is 

-9 . 9320240E-02 

The 

coef f 

of 

z 

is 

-1.295335 

The 

constant 

d 

is 

0.4731036 
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The input file was "STAMED52 .COD " 

The output file is ■ STAMED52 . COE " 

The coeff of x-squared is -4 .7458898E-02 
The coeff of y-squared is 0.3875865 

The coeff of z-squared is 0.6523174 

The coeff of yz is -0.7277786 

The coeff of zx is 0.3749601 

The coeff of xy is -0.4168406 

The coeff of x is -0.5613320 

The coeff of y is 0.9560621 

The coeff of z is -1.280479 

The constant d is 0.7334089 
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Best fit plot for the averaged sphere image. 


48 




33111111111111111 131333333333333333355151 
33311111111111111 131333333333333333333151 
15111111111111111111113333333333333333331 
15111111111111111111113133333333333333331 
31331111111111311 13333313133333333333353 1 
31311111111111111 133333 1 3133333333333553 1 
15111111113111111111113133333333553555551 
151111113131111111111 13333333333333335551 
3511111111111 1313133335555353333333555555 
351111111111113111331 35555553555355555555 
5531111111311111111 155553533335555555555 1 
553511 1 1 1 1311311 1 155555555353355555555551 
35313311111111113155151155555555555555551 
353533131 111111 15555555555555555555555551 
555555 11311 15155555533335555555555555555 1 
55555555551151555555333355555555555555551 
55353555555555555555555555555555355555551 
55555555555555555555555555555555555555551 
55555555355355555555555555555555555555551 
55555555355555555555555555555555555555551 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555551 
55555555555555555555555555555555555555551 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 
55555555555555555555555555555555555555555 


Best fit plot for the cylinder raw image. 
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Best fit plot for the averaged cylinder image. 
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SURFACE CHARACTERIZATION USING 3-D DISCRIMINANT APPROACH FOR SIMULATED DATA 


CC 

EFFICIENTS 

OF THE SI 

MULATED 

OBJECTS 




7 

A, OOEFF . OPX 

1 

0 

0 

1 

0 

0 

1 

H 


4 

1 

0 

3 

0 

0 

0 

0 

C. COEFF. OF Z L 

9 

20 

0 

2 

6 

1 

3 

2 

P, COEFF. OP YZ 

-6 

-4J3 

1J 

•1 


3 

1 


O. OOEIT. OP XZ 

3 

-2-5 

1 

0 

1 

EH 

•1 


H, COEW. OFXY 

mm 

OJ 

OJ 

1 

OJ 

l 

0 


P. COEFF. OF X 

i 

OJ 

*1 

2 

am 

2 

0 

0 

Q. COEFF. OP Y 

7 

0 

0 

3 

3 

3 

1 

■ 

R, COEFF . OF Z 

0 

0 

3 

0 

0 

0 

3 

H 

D. CONSTANT 

10 

0 

0 

8 

0 

12 

9 

19 

OBJECT IS 

(» 

( 2 ) 

(3) 

(4) 

(3) 

(« 

(7) 

»> 


Object 

(1) refers to PARABOLIC CYLINDER 

(2) refers to HYPERBOLIC PARABOLOID 

(3) refers to HYPERBOLOID OF ONE SHEET 

(4) refers to ELLIPSOID 

(5) refers to HYPERBOLIC CYLINDER 

(6) refers to REAL QUADRIC CONE 

(7) refers to HYPERBOLOID OF TWO SHEETS 

(8) refers to ELLIPTIC PARABOLOID 
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c************** SAMPLE DATA OF 3-D DISCRIMINANT PROGRAM 


Coeff. of x A 2 (A): 
A = 


1 


Coeff. of y A 2 (B): 
B = 


4 


Coeff. of z A 2 (Q: 
C = 


9 


Coeff. of yz (F): 
F = 

-6 


Coeff. of xz (G): 
G = 


3 


Coeff. of xy (H): 
H = 


-2 


Coeff. of x (P): 
P = 


1 


Coeff. of y (Q): 
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Q * 


7 


Coeff. of z (R): 
R = 


0 


Constant of prop. (D): 
D = 

10 


e = 


1 -2 3 

-2 4-6 

3-6 9 


1-231 
-2 4-6 7 

3-690 
1 7 0 10 


dt_e = 


0 


dt_EE * 
0 
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KJ C* 

- 0.0000 

0.0000 

14.0000 

rho_3 = 

1 

rho_4 = 

3 

s_d„ EE a 
0 

si = 

~1 

S2 as 

1 

S3 a 

1 

flag a 
0 

^ S,g ” ° f lhe ch - "*><* «« not the same 
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The rank of EE is : 1.0000 

The rank of e is : 3.0000 

The sign of the determinant of EE is : 0.0000 

The characterstics roots have the same sign? : 
The object is a PARABOLIC CYLINDER 


0.0000 



APPENDIX C 



C***** implementation of the 3-D discriminant approach 
diary on 

inputfCoeff. of x A 2 (A): ’); 

A=ans 

input(’Coeff. of y A 2 (B): ’); 

B=ans 

input(’Coeff. of z A 2 (C): ’); 

C=ans 

input(’Coeff. of yz (F): ’); 

F=ans 

input(’Coeff. of xz (G): ’); 

G=ans 

inputfCoeff. of xy (H): ’); 

H=ans 

input(’Coeff. of x (P): ’)*. 

P=ans 

input(’Coeff. of y (Q): ’); 

Q=ans 

input(’Coeff. of z (R): ’); 

R=ans 

input(’ Constant of prop. (D): ’); 

D=ans 

F=F/2; 

G=G/2; 

H=H/2; 

P=P/2; 

Q=Q/2; 

R=R/2; 

e=[A H G 
H B F 
G F C] 

EE=( A H G P 
H B F Q 
G F C R 
P Q R D ] 

dt_e=det(e) 

dt_EE=det(EE) 

K_K=eig(e) 

rho_3=rank(e) 

rho_4=rank(EE) 

s_d_EE=sign(dt_EE) 

sl=sign(K_K(l)> 

s2=sign(K_K(2)) 
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s3=sign(K_K(3)) 

flag=0 

if si = s2; 

flag=flag+l 

end; 

if si = s3; 
flag=flag+l 
end; 

if flag = 2; 
an_w=l; 

fprintf(’\n\n The sign of the ch. roots are the same \n’) 

else; 

an_w=0; 

end; 

fprintf(’\n\n The sign of the ch. roots are not the same \n’> 


fprintf(’\n\n The rank of EE is : %9.4f \n rho_3 ) 

fprintffVNn The rank of e is : %9.4f \n \ rho_4 ) 

fjjrintfCVnNn The sign of the determinant of EE is : %9.4f \n\s_d_EE ) 

fprintf(’\n\n The characterstics roots have the same sign? : %9.4f \n\ an_w) 

if rho_3=3 
if rho_4=4 
if s_d_EE=- 1 
if an_w==l 

fprintf(’\n\n The object is an ELLIPSOID Nn\n’) 

end 

end 

end 

end 

if rho_3=3 
if tho_4 = 4 
if s_d_EE=l 
if an_w=0 

fprintf(’\n\n The object is a HYPERBOLOID OF ONE SHEET \n\n’) 

end 

end 

end 

end 

if rho_3=3 
if rho_ 4 — 4 
if s_d_EE=- 1 
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if an_w=0 

fprintf(’\n\n The object is a HYPERBOLOID OF TWO SHEETS \n’) 

end 

end 

end 

end 

if rho_3=3 
if rho_4=3 
if an_w=0 

fprintfCNnNn The object is a REAL QUADRIC CONE \n’) 

end 

end 

end 

if rho_3=2 
if rho_4 — 4 
if s_d_EE==-l 
if an_w==0 

fprintf(’\n\n The object is an ELLIPTIC PARABOLOID \n’) 

end 

end 

end 

end 

if iho__3— 2 
if fho _ 4 — 4 
if s_d_EE=l 
if an_w=0 

fprintf(’\n\n The object is a HYPERBOLIC PARABOLOID V) 

end 

end 

end 

end 

if rho_3— 2 
if rho_4=3 
if an_w==l 

fprintf(’V\n The object is an ELLIPTIC CYLINDER Nn’) 

cad 

end 

end 

if rho_3— 2 
if rho_4=3 
if an_w=M) 

fprintfC’VnNn The object is a HYPERBOLIC CYLINDER \n’) 

end 

end 

end 
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if rho_3=l 
if rfio_4=3 

fprintfCVnNn The object is a PARABOLIC CYLINDER V) 

end 

end 

diary off 
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c***** Implementation of the 3-D discriminant ap poach 
diary on 

inputfCoeff. of x A 2 (A): ’); 

A=ans 

input(’Coeff. of y A 2 (B): ’); 

B=ans 

input(’Coeff. of z A 2 (Q: ’); 

C=ans 

inputfCoeff. of yz (F): ’); 

F=ans 

input(’Coeff. of xz (G): ’); 

G=ans 

input(’Coeff. of xy (H): ’); 

H=ans 

input(’Coeff. of x (P): ’); 

P=ans 

input(’Coeff. of y (Q): ’); 

Q=ans 

input(’Coeff. of z (R): ’); 

R=ans 

input(’Constant of prop. (D): ’); 

D=ans 

F=F/2; 

G=G/2; 

H=H/2; 

P=P/2; 

0 = 0 / 2 ; 

R=R/2; 

e=[A H G 
H B F 
G F C] 

EE=[ A H G P 
H B F Q 
G F C R 
P Q R D ] 

dt_e=det(e) 

dt_EE=det(EE) 

K_K=eig(e) 

rho_3=rank(e) 

rho_4=rank(EE) 

s_d_EE=sign(dt_EE) 

sl=sign(K_K(l)) 

s2=sign(K_K(2)) 
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s3=sign(K_K(3)) 

flag=0 

if si = s2; 

flag=flag+l 

end; 

if si = s3; 
flag=flag+l 
end; 

if flag = 2; 
an_w=l; 

fprintf(’\n\n The sign of the ch. roots are the same \n’) 

else; 

an_w=0; 

end; 

fprintf(’\n\n The sign of the ch. roots are not the same \n’> 


fprintf(’\ri\n The rank of EE is : %9.4f \n rho_3 ) 

fprintf(’\n\n The rank of e is : %9.4f \n \ rho_4 ) 

fprintf(’\n\n The sign of the determinant of EE is : %9.4f \n\s_d_EE ) 

fprintf(V\n The characterstics roots have the same sign? : %9.4f \n\ an_w) 

if rfio_3=3 
if iho_4=4 
if s_d_EE=- 1 
if an_w=l 

fprintfCWt The object is an ELLIPSOID \n\n’) 

end 

end 

end 

end 

if rho_3=3 
if rho_4==4 
if s_d_EE=l 
if an__w==0 

fprintf(’\ri\n The object is a HYPERBOLOID OF ONE SHEET \n\n’) 

end 

end 

end 

end 

if rho_3=3 
if rho_4=4 
if s_d_EE=-l 


61 



if an_w=0 

fprintf(’\n\n The object is a HYPERBOLOID OF TWO SHEETS \n’> 

end 

end 

end 

end 

if rho_3— 3 
if rho_4=3 
if an_w=0 

fprintf(’\n\n The object is a REAL QUADRIC CONE \n’) 

end 

end 

end 

if rho_3=2 
if rho_4==4 
if s_d_EE=-l 
if an_w=0 

fprintf(’\n\n The object is an ELLIPTIC PARABOLOID \n’) 

end 

end 

end 

end 

if rho_3==2 
if rho_4=4 
if s_d_EE=l 
if an_w=0 

fprintf(’\n\n The object is a HYPERBOLIC PARABOLOID \n’) 

end 

end 

end 

end 

if rho_3— 2 
if rho_4=3 
if an_w=l 

fprintfCVnNn The object is an ELLIPTIC CYLINDER \n’) 

end 

end 

end 

if rho_3=2 
if rho_4==3 
if an_w=0 

fprintf(’Nn\n The object is a HYPERBOLIC CYLINDER \n’) 

end 

end 

end 
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if rho_3=l 
if rho_4=3 

fprintf( , NnNn The object is a PARABOLIC CYLINDER \n’) 

end 

end 

diary off 
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c**** PROGRAM SURFACE ALIGNMENT 

c**** This program is used to ged rid of the product terms 

c**** from the quadratic representation of any 3D surface. 

c**** The new coefficents generated consisits of the square terms 

c**** and the x,y ,z, and the constant term. 

c**** Last compiled on 14feb. 1991. 

REAL AA,BB,CC,DD,FF,GG,HH,PP,QQ,RR,D,Test_f,Test_g,test_h 
REAL A(50,50),B(50,50),C(50,50),F(50,50) 

REAL G(50,50),H(50,50),ALPHA(100)3ETA(100) 

REAL RESULT(200,200),P(50,50),Q(50,50) 1 R(50,50) 

REAL AAA,BBB,CCC,DDD,EEE,FFF,GGG3HHjn,ROT(3,3) 

REAL DEL1,DEL2,DEL3A_A3_B,C_C3_F,G_G3_H,GAMMA(100) 
REAL W,WV,VWV,VWW,THRESHLDJNITMIN,ABSA,ABSB,ABSC 
REAL A_AA,B_BB,C_CC,D_DD,P_PP,Q_QQ3_RR 
REAL ABSF,ABSG,ABSH,ABSP,ABSQ,ABSR,RRR(50),alptot,bettot 
REAL gamtot 
INTEGER NJtfJU 

C F(X,Y,Z)=Ax**2+By**2+Cz**2+2Fyz+2Gxz+2Hxy+2Px+2Qy+2Rz+D 
C =0 

C PARAMETER (THRESHLD= 0.00000000000000001) 

OPEN(UNIT=l,FILE=’CONVERGENCE.DAT’,STATUS=’NEW’) 

TYPE*, ’ENTER VALUE FOR THRESHLD:’ 

ACCEPT* .THRESHLD 

Type*,’Enter coef. of x ** 2 
Accept* AA 

Type*,’Enter coef. of y ** 2 
Accept*3B 

Type*,’Enter coef. of z ** 2 
Accept*, CC 

Type*, ’Enter coef. of yz 
Accept* 3F 

Type*,’Enter coef. of xz :’ 

Accept* ,GG 

Type*, ’Enter coef. of xy :’ 

Accept* 3H 

Type*, ’Enter coef. of x 

Accept*3P 

Type*, ’Enter coef. of y :’ 

Accept*,QQ 

Type*, ’Enter coef. of z :’ 

Accept* 3R 

Type*, ’Enter constant of prop. :’ 

Accept* JD 
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A(1,1)=AA 

B(1,1)=BB 

C(1,1)=CC 

F(1,1)=FF 

G(1,1)=GG 

H(1,1)=HH 

P(1,1)=PP 

Q(U)=QQ 

R(1,1)=RR 

ABSA=ABS(AA) 

ABSB=ABS(BB) 

ABSC=ABS(CC) 

ABSF=ABS(FF) 

ABSG=ABS(GG) 

ABSH=ABS(HH) 

ABSP=ABS(PP) 

ABSQ=ABS(QQ) 

ABSR=ABS(RR) 

RRR(1)=ABSA 

RRR(2)=ABSB 

RRR(3)=ABSC 

RRR(4)=ABSF 

RRR(5)=ABSG 

RRR(6)=ABSH 

RRR(7)=ABSP 

RRR(8)=ABSQ 

RRR(9)=ABSR 

DO 3980 1=1,9 

IF (RRR(I).EQ.0)THEN 

RRR(I)= 10000 

ENDIF 

3980 CONTINUE 

INITMIN=AMIN1(RRR(1),RRR(2),RRR(3),RRR(4),RRR(5),RRR(6),RRR(7) 
+ JRRR(8),RRR(9)) 

WRITE(*,*)INrrMIN 

IF (ABS(INITMIN).LT.1.0)THEN 

A(i,i)=A(i,iyiNrrMiN 

B(1,1)=B(1,1)/INITMIN 

C(1,1)=C(1,1)/INITMIN 

F(1,1)=F(1,1)/INITMIN 

G(l,l)=G(l,iyiNITMIN 

H(1,1)=H(1,1)/INITMIN 

P(1,1)=P(1,1)/INITMIN 

Q(l,l)=Q(U)/INrrMIN 


65 



Q(l t l)=Q(UyiNITMIN 

DD_D=D/INITMIN 

ELSE 

GOTO 3405 
ENDIF 

3405 A(1,1)=AA 
B(1,1)=BB 
C(1,1)=CC 
F(1,1)=FF 
G(1,1)=GG 
H(1,1)=HH 
P(1,1)=PP 
Q (U)=QQ 
R(1,1)=RR 


345 if (b(l,l).eq.a(l,l)) then 

goto 1167 
else 

c goto 57 

c endif 

c else 

goto 57 
endif 

57 alpha(l)=(0.5*ATAND((H(U)/(B(U)-A(U))») 

A(U)=A(l,l)*COSD(ALPHA(l))*COSD(ALPHA(l))+B(l,l)* 

+ SIND(ALPH A(1))*SIND( ALPHA(1 ))- H(1,1)*SIND(ALPHA(1))* 
+ COSD( ALPHA( 1 » 

B(1 , 1 )=B( 1 , 1 )*COSD( ALPHA( 1 ))*COSD(ALPHA( 1 ))+A( 1,1)* 

+ SIND(ALPHA(1))*SIND(ALPHA(1))+H(1,1)*SIND(ALPHA(1))* 
+ COSD(ALPHAG)) 

CCUMXU) 

F( 1 , 1 )=G( 1 , 1)*SIND( ALPHA( 1 ))+F(l , 1 ) *COSD( ALPHA( 1 )) 
G(l,l)=G(l,l)*COSD(ALPHA(l))-F(l,l)*SIND(ALPHA(l)) 
H(1,1)=0 

P(l,l)=P(U)*COSD(ALPHA(l))-Q(l,l)*SIND(ALPHA(l)) 

Q(U)=Q(U)*COSD(ALPHA(l))+P(l,l)*SIND(ALPHA(l)) 

R(1,1)=R(1,1) 


IF ( ABS(F(1 , 1)).LT.THRESHLD)THEN 
GOTO 1005 
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ELSE 

GOTO 1167 
END IF 

1005 IF (ABS(G( 1 , 1)) LT.THRESHLD)THEN 
GOTO 1812 
ELSE 

GOTO 1167 
END IF 

1167 IF (C(1,1).EQ.B(1,1))THEN 
GOTO 1169 
ELSE 

GOTO 1200 
ENDIF 

1200 BETA(1)=(0.5*ATAND((F(1,1)/(C(1,1)-B(1,1))))) 

A(U)=A(1,1) 

B(l,2)=B(l,l)*COSD(BETA(l))*COSD(BETA(l))+C(l,l)* 

+ SIND(BETA(l))*SIND(BETA(l))-F(l,l)*SIND(BETA(l))*COSD(BETA(l)) 


C(l,2)=C(l,l)*COSD(BETA(l))*COSD(BETA(l))+B(l,l)* 

+ SIND(BETA(l))*SIND(BETA(l»+F(l,l)*SIND(BETA(l))*COSD(BETA(l)) 
F(U)=0 

G(U)=G(l,l)*COSD(BETA(l)) 

H(U)=-G(1,1)*SIND(BETA(1)) 

P(U)=P(1,1) 

Q(U)=Q(U)*COSD(BETA(l))-R(U)*SIND(BETA(l)) 

R( 1 ,2)=R( 1 , 1 )*COSD(BET A( 1 ))+Q( 1 , 1 )*SIND(BETA( 1 )) 

IF (ABS(H(U)).LT.THRESHLD)THEN 

GOTO 1007 

ELSE 

GOTO 1169 
ENDIF 

1007 IF (ABS(G(U)).LT.THRESHLD)THEN 
GOTO 1812 
ELSE 

GOTO 1169 
ENDIF 


1169 IF (C(U).EQ.A(U))THEN 
GOTO 67 
ELSE 

GOTO 1235 
ENDIF 
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1235 GAMMA(1)=(0.5* ATAND((G(1 ,2)/(C(U)- A(1 ,2))))) 

A( 1 ,3)= A( 1 ,2)*C0SD(GAMMA( 1 ))*COS D(G AMMA( 1 ))+€( U) * 

+ SIND(GAMMA(1))*SIND(GAMMA(1))-G(U)*SIND(GAMMA(1)) 
+ *COSD(GAMMA(l» 

B(1,3)=B(U) 

C(1,3)=C(U)*C0SD(GAMMA(1))*C0SD(GAMMA(1))+A(U)* 

+ SIND(GAMMA(1))*SIND(GAMMA(1))4G(U)*SIND(GAMMA(1)) 
+ *C0SD(GAMMA(1)) 

F(1,3)=H(U)*SIND(GAMMA(1)) 

G(l,3)=0 

H( 1 ,3)=H( 1 ,2)*C0SD(G AMMA( 1 )) 

P(l,3)=P(l,2)*COSD(GAMMA(l))-R(l,2)*SIND(GAMMA(l)) 

QCUHXU) 

R(1,3)=R(U)*C0SD(GAMMA(1))+P(U)*SIND(GAMMA(1)) 


IF (ABS(F(1,3)).LT.THRESHLD)THEN 

GOTO 1009 

ELSE 

GOTO 67 

END IF 

1009 IF (ABS(H(1,3)).LT.THRESHLD)THEN 
GOTO 1812 
ELSE 
GOTO 67 
ENDIF 


67 DO 10 1=2,100 

C DO 20 J=1 

71 if((b(i- 1 ,3).eq.a(i- 1 ,3)))then 

goto 167 
C else 

c if(h(i,3).eq.0)then 

c goto 67 

c else 

c goto 67 

c endif 

else 

goto 177 
endif 
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177 alpha(D=(0.5*ATAND((H(I- 1 ,3)/(B(I- 1 ,3)-A(I- 1 3))))) 

Aa,l)=Aa- 1 ,3)*C0SD(ALPHA(I))*C0SD(ALPHA(I))+(B(I- 1 ,3))* 
+ SIND(ALPHA(I))*SIND(ALPHA(I))- Ha-1,3)*SIND(ALPHA(I))* 
+ COSD(ALPHA(I)) 

B(I,l)=B(I-l,3)*COSD(ALPHA(I))*COSD(ALPHA(I))+A(I-l,3)* 

+ SIND( ALPHA (D)*SIND(ALPHA(I))+H(I- 1 ,3)*SIND(ALPHA(I))* 
+ COS D( ALPHA(D) 
c(i,D=ca-i.3) 

Fa, l)=Fa-l,3)*COSD(ALPHA(D) 

Ga,l)=-Fa- l,3)*SIND(ALPHAa)) 

Ha,i)=o 

pa, l)=Pa- 1 ,3)*COSD( ALPHA0))-Qa- 1 ,3)*SIND( ALPHAO)) 
Qa,l)=Qa- 1 ,3)*COSD(ALPHA0))+Pa- 1 ,3)*SIND( ALPHAO)) 
RO,l)=Ra-l,3) 


IF (ABSO J a,l)) LT.THRESHLD)THEN 

GOTO 1011 

ELSE 

GOTO 167 

ENDIF 

1011 IF (ABS(Gai))-LT.THRESHLD)THEN 
N=I 

GOTO 666 
ELSE 
GOTO 167 
ENDIF 

167 if «c(i, 1 ).eq.b(i, l)))then 

goto 69 
else 

c goto 59 

c endif 

c else 

goto 59 
endif 

59 BETA(I)=(0.5*ATAND((Fa,l)/(Ca.l)-Bai))))) 

A0,2)=A(U) 

Baa)=Ba,l)*COSD(BETA0))*COSD(BETAa»K:a,l)* 

+ SIND(BETAa))*SIND(BETAa))-Fa,l)*SIND(BETAa))*COSD(BETAa)) 

Ca^)=Ca,l)*COSD(BETA0))*COSD(BETA0)>t-Ba,l)* 

+ SIND(BETAO»*SIND(BETAa))+Fa, l)*SIND(BETAa»*COSD(BETAa)) 
F0,2)=0 

G0^)=Ga l)*COSD(BETA0)) 

Ha^)=-G0, 1)*SIND(BETA0)) 

P0,2)=Pa,l) 

Q0^)=CKU)*COSD(BETAa))-Ra,l)*SIND(BETAa)) 
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R(L2)=R(I,1)*C0SD(BETA(I))+Q(M)*SIND(BETA(I)) 


IF (ABS(G(I,2)).LT.THRESHLD)THEN 

GOTO 1013 

ELSE 

GOTO 69 

ENDIF 

1013 IF (ABS(H(I,2)).LT.THRESHLD)THEN 
N=I 

GOTO 666 
ELSE 
GOTO 69 
END IF 

c69 if(g(i,2).eq.0)then 
69 tf((c(U).eq.a(U)))then 

goto 10 
else 

c goto 63 

c endif 

c else 

goto 63 
endif 

63 GAMMA (I)=(0.5*ATAND((G(I,2)/(C(L2)-A(I t 2))))) 

A(I,3)=A(I,2)*COSD(GAMMA(I))*COSD(GAMMA(D)+C(I,2)* 

+ SIND(GAMMA(I))*SIND(GAMMA(I))-G(I,2)*SIND(GAMMA(I)) 
+ *COSD(GAMMA(I)) 

B(L3)=B(U) 

C(I,3)=C(I.2)*COSD(GAMMA(I))*COSD(GAMMA(I))+A(I,2)* 

+ SIND(GAMMA(I))*SIND(GAMMA(I))-K}(I,2)*SIND(GAMMA(I)) 
+ *COSD(GAMMA(I» 

F(L3)=H(I,2)*SIND(GAMMA(I)) 

G(I,3)=0 

H(I,3)=H(I,2)*COSD(GAMMA(I)) 

Pa3)=P(L2)*COSD(GAMMA(I))-R(I,2)*SIND(GAMMA(I)) 

qcumxu) 

R(L3)=R(I,2)*COSD(GAMMA(I))+P(I,2)*SIND(GAMMA(I)) 

IF (ABS(F(I,3)).LT.THRESHLD)THEN 

GOTO 1015 

ELSE 

GOTO 10 

ENDIF 

1015 IF (ABS(H(L3)).LT.THRESHLD)THEN 
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N=I 

GOTO 666 
ELSE 
GOTO 10 
END IF 


20 CONTINUE 
10 CONTINUE 
1812 N=1 
666 WRTTE(**)N 
WRITER, 123) 

123 FORMAT^X,’**************************************************** 

+ ************>) 

WRITE(*,*)(’THE NUMBER OF ITERATIONS COMPLETED IS:\N) 
M=N*3 

DO 1000 1=1, N 
DO 1001 J=l,3 
RESULT(3*(I- 1)+J, 1)=A(IJ) 

RESULT(3*(I- 1)+L2)=B(IJ) 

RESULT(3*(I-1)+L3)=C(IJ) 

RESULT(3*(I- 1)+J,4)=F(IJ) 

RESULT(3*(I-0+J.5)=G(U) 

RESULT(3*(I- 1)+J,6)=H(I J) 

RESULT(3*(I- 1)+J,7)=P(U) 

RESULT(3*(I-1)+L8)=Q(IJ) 

RESULT(3*(I-1)+L9)=R(I4) 

1001 CONTINUE 
1000 CONTINUE 

WRTTE(1,*X’THE NUMBER OF ITERATIONS COMPLETED IS :\N) 
WRTTE(1,123) 


WRrrE(l,*X’COEFF. OF X SQUARE TERM IS : \ AA) 
198 WRrrE(l,*X’COEFF. OF Y SQUARE TERM IS : BB) 

298 WRTTE(l,*X’COEFF. OF Z SQUARE TERM IS : CC) 

398 WRrrE(l,*)0COEFF. OF YZ SQUARE TERM IS : FF) 

498 WRITE( 1 ,*)(’ COEFF. OF XZ SQUARE TERM IS : ’ ,GG) 
598 WRrrE(l,*)(’COEFF. OF XY SQUARE TERM IS : ’ 3H) 
WRTTECl^X’COEFF. OF X TERM IS : \ PP) 

WRITE(1 ,*X’COEFF. OF Y TERM IS : QQ) 

WRTTE(1 ,♦)(’ COEFF. OF Z TERM IS : RR) 

WRrrE(l,*)(’CONSTANT OF PROP. IS : \ D) 
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write(l,123) 

write(l,123) 

write(l,123) 

DO 2000 I=1,M 

WRITE(1,*)(RESULT(I,J)J=1,9) 
2000 CONTINUE 

A_AA=RESULT (M-2, 1 ) 

B_BB=RESULT(M-2,2) 

C_CC=RESULT(M-2,3) 

PJPP=RESULT(M-2,7) 

Q_QQ=RESULT(M-2,8) 

R_RR=RESULT(M-2,9) 

D_DD=D 


do 30001 i=l,3 
write(l,123) 

30001 continue 

WRrrE(l,*)(’THE NEW COEFF. OF X SQUARE TERM IS : \ A_AA) 
WRITE(1,*)(’THE NEW COEFF. OF Y SQUARE TERM IS : \ B_BB) 
WRrTE(l,*)(’THE NEW COEFF. OF Z SQUARE TERM IS : \ C_CQ 
WRTTE(1,*)(THE NEW COEFF. OF X TERM IS : \ P_PP) 
WRITE(1,*)(’THE NEW COEFF. OF Y TERM IS : \ Q_QQ) 
WRrrE(l,*)(’THE NEW COEFF. OF Z TERM IS : \ R_RR) 
WRrrE(l,*)(’THE NEW CONSTANT OF PROP. IS : \D_DD) 

do 3001 i=l,3 
write(l,123) 

3001 continue 

write(l,1278) 

1278 formatCex/A’^/B’^’C’^’F’^’G’^/H’^’P’, 

+ 9x,’Q’,9x,’R’) 
write(l,1897) 

1897 format(5x,’ 

+ ’) 


DO 2001 1=1 M 

WRITE(1,1234)(RESULT(IJ)J= 1 ,9) 
2001 CONTINUE 
1234 format(9F10.5) 

DO 3000 1=1,5 
WRTTE(1,123) 

3000 CONTINUE 
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/ 


write(l,1908) 

1908 format(6x, ’Alpha’, 9x,’Beta’,9x,’Gamma’) 
write(l,1897) 

DO 4000 I=1,N 

WRITE(1,*)ALPHA(1)3ETA(I),GAMMA(I) 

4000 CONTINUE 

alptot=alpha(l)+alpha(2)+alpha(3) 

bettot=beta(l)+beta(2)+beta(3) 

gamtot=gamma( 1 )+gamma(2)+gamma(3) 

write(l,123) 

write(l,123) 

write(l,1998) 

1998 format(6x,’ALPTOT’,9x,’BETTOT’,9x,’GAMTOT’) 
write(l ,*)alptot,bettot,gamtot 
write(l,123) 

c***** To evaluate coeff. of yz, xz, and xy once alpha, beta 
c***** and gamma are evaluated. 

write(*,*)alpha( l),beta( 1 ),gamma( 1 ) 

A AA=BB *cosd(alpha( 1 )) *cosd(alpha( 1 ))+(AA *sind(alpha( 1 )) 

+ *sind(alpha(l)))+((HH/2) # sind(2*alpha(l)))-CC 
BBB=gg*sind(alpha(l))+(ff*cosd(alpha(l))) 
CCC=((aa-bb))*sind(2*alpha(l))+(hh*cosd(2*alpha(l))) 
DDD=gg*cosd(alpha(l)Hff*sind(alpha(l») 
EEE=aa*(cosd(alpha(l))*cosd(alpha(l»(sind(alpha(l» 

+ *sind(alpha(l))*sind(beta( l))*sind(beta(l)») 

FFF=bb*(sind(alpha(l))*sind(alpha(l))-(cosd(alpha(l» 

+ *cosd(alpha( 1 ))*sind(beta( 1 ))*sind(beta( 1 )))) 

GGG=cc*cosd(beta( 1 ))*cosd(beta( 1)) 

HHH=(gg/2)*sind(alpha( 1 ))*sind(2*beta( 1 ))+((ffZ2)*cosd(alpha(l))* 

+ sind(2*beta(l))) 

HI=(hh/2)*sind(2*alpha( 1 ))*(1 +sind(beta( 1 ))* 

+ sind(beta(l») 

Test_F=(AAA*sind(2*beta(l))+BBB*cosd(2*beta(l)))* 

+ cosd(gamma( l))+(CCC*cos(beta( 1 ))-DDD*sind(beta(l))) 

+ *and(gamma(l)) 

Test_G=(EEE+FFF-GGG-HHH-m)*SIND(2*GAMMA(l)) + 

+ (CCC*SIND(BETA(l))+COSD(BETA(l))*DDD)*COSD(2*GAMMA(l)) 

TEST_H=(CCC*COSD(BET A( 1 )>DDD* SIND<BETA( l»)*COSD(G AMMA(1 )) 

+ (AAA*SIND(2*BETA(l»+BBB*COSD(2*BETA(l)))*SIND(GAMMA(l» 

c write(l,123) 
write(l,123) 
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write(l,124) 


c 
c 

c write( 1 ,*)test_f,test_g,test_h 
write(l,123) 
write(l,123) 
write(l,124) 

124 format(5x,’THE ROTATION MATRIX IS : ’) 
c**** To evaluate the rotation matrix 

rot(l,l)=cosd(alpha(l))*cosd(gamma(l))-(sind(alpha(l»* 
+ sind(beta( 1 )) * sind(gamma( 1 ))) 

rot(l,2)=-sind(alpha(l))*cosd(gamma(l»-(cosd(alpha(l»* 
+ sind(beta(l))*sind(gamma(l))) 

rot( 1 ,3)=-sind(gamma( 1 )) *cosd(beta( 1 )) 

rot(2, 1 )=sind(alpha( 1 )) *cosd(be ta( 1 )) 
rot(2,2)=cosd(beta( 1 ))*cosd(alpha( 1 )) 
rot(2,3)=-sind(beta(l )) 

rot(3,l)=cosd(alpha(l))*sind(gamma(l»+(sind(alpha(l))* 
+ sind(beta(l))*cosd(gamma(l))) 

rot(3,2)=cosd(alpha(l))*cosd(gamma(l))*sind(beta(l)) 

+ -(sind(alpha(l))*sind(gamma(l)» 
rot(3,3)=cosd(gamma( l))*cosd(beta(l» 

DO 989 1=1,3 

WRITE(l,*)(ROT(I,J)J=l,3) 

989 CONTINUE 
stop 
end 
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THE COEFFICIENTS CONSIDERED ARE OF THE AVERAGED CYLINDER. 
********************************************************************* 

THE NUMBER OF ITERATIONS COMPLETED IS : 4 

********************************************************************* 

COEFF. OF X SQUARE TERM IS : 0.8338000 
COEFF. OF Y SQUARE TERM IS : 4.1 100001E-03 
COEFF. OF Z SQUARE TERM IS : 5.9000000E-02 
COEFF. OF YZ SQUARE TERM IS : -1.0300000E-03 
COEFF. OF XZ SQUARE TERM IS : -0.6368000 
COEFF. OF XY SQUARE TERM IS : 0.4437000 
COEFF. OF X TERM IS : 1.1410000E-02 

COEFF. OF Y TERM IS : -0.1894000 
COEFF. OF Z TERM IS : 0.1932000 
CONSTANT OF PROP. IS : -0.1341000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 

0.8893950 -4.8199981E-02 5.9000000E-02 0.1537948 -0.5803153 

0.0000000E4O0 -3.497 1744E-02 -0.1752182 0.1932000 

0.8893950 -8.8334545E-02 9.9134572E-02 O.OOOOOOOE+OO -0.5144598 

0.2685087 -3.497 1744E-02 -0.2447266 9.0202741E-02 

0.9657465 -8.8334545E-02 2.2782966E-02 7.6404691E-02 O.OOOOOOOE+OO 

0.2574087 -5.9193403E-02 -0.2447266 7.6522537E-02 

0.9812339 -0.1038219 2.2782966E-02 7.5857453E-02 9.1281543E-03 

0.0000000E+O0 -8.8007204E-02 -0.2359019 7.6522537E-02 

0.9812339 -0.1143151 3.3276103E-02 O.OOOOOOOE+OO 8.7976847E-03 

-2.4339 152E-03 -8.8007204E-02 -0.2477653 1.0851689E-02 

0.9812542 -0.1143151 3.3255693E-02 1.1293819E-05 O.OOOOOOOE+OO 

-2.4338888E-03 -8.7955907E-02 -0.2477653 1.1259941E-02 

0.9812556 -0.1143164 3.3255693E-02 1.1293812E-05 -1.2545008E-08 

O.OOOOOOOE-KX) -8.7680638E-02 -0.2478629 1.1259941E-02 

0.9812556 -0.1143164 3.3255693E-02 0.0000000E-KX) -1.2545008E-08 

4.8003978E-13 -8.7680638E-02 -0.2478633 1.1250457E-02 

0.9812556 -0.1143164 3.3255693E-02 3.1762148E-21 O.OOOOOOOE+OO 

4.8003978E-13 -8.7680638E-02 -0.2478633 1.1250456E-02 

0.9812556 -0.1143164 3.3255693E-02 3.1762148E-21 6.9585087E-34 

O.OOOOOOOE-KX) -8.7680638E-02 -0.2478633 1.1 250456E-02 

0.0000000E+00 0.0000000E-KX) O.OOOOOOOE-KX) O.OOOOOOOE+OO O.OOOOOOOE+OO 

O.OOOOOOOE+OO 0.0000000E-KX) O.OOOOOOOE+OO O.OOOOOOOE+OO 

O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 

O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
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**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 


THE NEW COEFF. OF X SQUARE TERM IS : 0.9812556 
THE NEW COEFF. OF Y SQUARE TERM IS : -0.1143164 
THE NEW COEFF. OF Z SQUARE TERM IS : 3.3255693E-02 
THE NEW COEFF. OF X TERM IS : -8.7680638E-02 
THE NEW COEFF. OF Y TERM IS : -0.2478633 
THE NEW COEFF. OF Z TERM IS : 1.1250456E-02 
THE NEW CONSTANT OF PROP. IS : -0.1341000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 

ABCFGHPQR 

0.88939 -0.04820 0.05900 0.15379 -0.58032 0.00000 -0.03497 -0.17522 0.19320 

0.88939 -0.08833 0.09913 0.00000 -0.51446 0.26851 -0.03497 -0.24473 0.09020 

0.96575 -0.08833 0.02278 0.07640 0.00000 0.25741 -0.05919 -0.24473 0.07652 

0.98123 -0.10382 0.02278 0.07586 0.00913 0.00000 -0.08801 -0.23590 0.07652 

0.98123 -0.11432 0.03328 0.00000 0.00880 -0.00243 -0.08801 -0.24777 0.01085 

0.98125 -0.11432 0.03326 0.00001 0.00000 -0.00243 -0.08796 -0.24777 0.01126 

0.98126 -0.11432 0.03326 0.00001 0.00000 0.00000 -0.08768 -0.24786 0.01126 

0.98126 -0.11432 0.03326 0.00000 0.00000 0.00000 -0.08768 -0.24786 0.01125 

0.98126 -0.11432 0.03326 0.00000 0.00000 0.00000 -0.08768 -0.24786 0.01125 

0.98126 -0.11432 0.03326 0.00000 0.00000 0.00000 -0.08768 -0.24786 0.01125 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 

Alpha Beta Gamma 

-14.06845 27.56112 16.53207 

-6.861581 15.46432 -0.2658640 

6.3643321E-02 2.1924460E-03 3.7910129E-07 
-1.2552463E-11 0.0000000E+O0 O.OOOOOOOE+OO 

**************************************************** ************ 

**************************************************** ************ 

ALPTOT BETTOT GAMTOT 
-20.86638 43.02764 16.26621 

**************************************************** ************ 
**************************************************** ************ 
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**************************************************** ************ 
**************************************************** ************ 

THE ROTATION MATRIX IS : 

0.9619108 0.1053204 - 0.2522604 

- 0.2154955 0.8599276 - 0.4626946 

0.1681945 0.4994319 0.8498697 
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THE COEFFICIENTS CONSIDERED ARE OF THE 3X3 FILTERED IMAGE OF THE 
AVERAGED CYLINDER. 

******************************************************************** 

THE NUMBER OF ITERATIONS COMPLETED IS : 4 

******************************************************************** 

COEFF. OF X SQUARE TERM IS : 0.6636000 
COEFF. OF Y SQUARE TERM IS : 2.0900000E-02 
COEFF. OF Z SQUARE TERM IS : -9.2399999E-02 
COEFF. OF YZ SQUARE TERM IS : -2.1900000E-02 
COEFF. OF XZ SQUARE TERM IS : -0.7604000 
COEFF. OF XY SQUARE TERM IS : 0.7727000 
COEFF. OF X TERM IS : 0.4242000 
COEFF. OF Y TERM IS : -0.2155600 
COEFF. OF Z TERM IS : 0.3749000 

CONSTANT OF PROP. IS : -0.2534000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 

0.8447758 -0.1276161 -9.2399999E-02 0.3030199 -0.5598051 

O.OOOOOOOE+OO 0.2925456 -0.3193743 0.3749000 

0.8447758 -0.2625377 4.2521641E-02 O.OOOOOOOE+OO -0.4180661 

0.3722935 0.2925456 -0.4878348 6.7580700E-02 

0.8959736 -0.2625377 -8.6761 117E-03 8.8566750E-02 O.OOOOOOOE+OO 

0.3616053 0.2680698 -0.4878348 0.1352357 

0.9235348 -0.2900990 -8.6761 117E-03 8.7555319E-02 1.3346768E-02 

O.OOOOOOOE+OO 0.1914931 -0.5226611 0.1352357 

0.9235348 -0.2967516 -2.0234007E-03 O.OOOOOOOE-fOO 1.3195274E-02 

-2.00523 18E-03 0.1914931 -0.5370466 5.5175520E-02 

0.9235819 -0.2967516 -2.070428 IE-03 1.4292762E-05 O.OOOOOOOE+OO 

-2.005 1808E-03 0.1918815 -0.5370466 5.3809207E-02 

0.9235827 -0.2967524 -2.070428 IE-03 1.4292757E-05 -1.1742504E-08 

0.0000000E+O0 0.1923226 -0.5368887 5.3809207E-02 

0.9235827 -0.2967524 -2.070428 IE-03 0.0000000E+O0 -1.1742504E-08 

2.8476928E-13 0.1923226 -0.5368900 5.3796187E-02 

0.9235827 -0.2967524 -2.070428 IE-03 1.8062406E-21 O.OOOOOOOE+OO 

2.8476928E-13 0.1923226 -0.5368900 5.3796187E-02 

0.9235827 -0.2967524 -2.070428 IE-03 1.8062406E-21 2.1074614E-34 

O.OOOOOOOE+OO 0.1923226 -0.5368900 5.3796 187E-02 

O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO 0.0000000E+00 O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
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**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 

THE NEW COEFF. OF X SQUARE TERM IS : 0.9235827 
THE NEW COEFF. OF Y SQUARE TERM IS : -0.2967524 
THE NEW COEFF. OF Z SQUARE TERM IS : -2.070428 IE-03 
THE NEW COEFF. OF X TERM IS : 0.1923226 
THE NEW COEFF. OF Y TERM IS : -0.5368900 
THE NEW COEFF. OF Z TERM IS : 5.3796 187E-02 

THE NEW CONSTANT OF PROP. IS : -0.2534000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 

A B C F G H P Q R 

0.84478 -0.12762 -0.09240 0.30302 -0.55981 0.00000 0.29255 -0.31937 0.37490 

0.84478 -0.26254 0.04252 0.00000 -0.41807 0.37229 0.29255 -0.48783 0.06758 

0.89597 -0.26254 -0.00868 0.08857 0.00000 0.36161 0.26807 -0.48783 0.13524 

0.92353 -0.29010 -0.00868 0.08756 0.01335 0.00000 0.19149 -0.52266 0.13524 

0.92353 -0.29675 -0.00202 0.00000 0.01320 -0.00201 0.19149 -0.53705 0.05518 

0.92358 -0.29675 -0.00207 0.00001 0.00000 -0.00201 0.19188 -0.53705 0.05381 

0.92358 -0.29675 -0.00207 0.00001 0.00000 0.00000 0.19232 -0.53689 0.05381 

0.92358 -0.29675 -0.00207 0.00000 0.00000 0.00000 0.19232 -0.53689 0.05380 

0.92358 -0.29675 -0.00207 0.00000 0.00000 0.00000 0.19232 -0.53689 0.05380 

0.92358 -0.29675 -0.00207 0.00000 0.00000 0.00000 0.19232 -0.53689 0.05380 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 

Alpha Beta Gamma 

-25.12386 41.68550 13.76233 

-8.667336 8.640893 -0.4083926 

4.7072496E-02 1.3894888E-03 3.6341686E-07 

-6.6850809E-12 0.0000000E+O0 0.0000000E+O0 

**************************************************** ************ 

**************************************************** ************ 

ALPTOT BETTOT GAMTOT 
-33.74412 50.32778 13.35394 

**************************************************** ************ 

**************************************************** ************ 
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**************************************************** ************ 
**************************************************** ************ 


THE ROTATION MATRIX IS : 
0.9465714 0.2691452 -0.1776615 

-0.3170765 0.6761528 -0.6650413 

-5.8866352E-02 0.6858414 0.7253664 
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THE COEFFICIENTS CONSIDERED ARE OF THE 5X5 FILTERED IMAGE OF THE RAW 
CYLINDER. 

*********************************************************************** 

THE NUMBER OF ITERATIONS COMPLETED IS : 4 

*********************************************************************** 

COEFF. OF X SQUARE TERM IS : 5.7200000E-02 
COEFF. OF Y SQUARE TERM IS : 0.5990000 
COEFF. OF Z SQUARE TERM IS : 0.4416000 
COEFF. OF YZ SQUARE TERM IS : -0.8076000 
COEFF. OF XZ SQUARE TERM IS : 0.4595000 
COEFF. OF XY SQUARE TERM IS : -0.1490000 
COEFF. OF X TERM IS : -0.5915000 
COEFF. OF Y TERM IS : 1.089800 

COEFF. OF Z TERM IS : -1.019800 

CONSTANT OF PROP. IS : 0.6643000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 

4.7142603E-02 0.6088773 0.4416000 -0.8618143 0.3400714 

O.OOOOOOOE+OO -0.4403836 1.138920 -1.019800 

4.7142603E-02 0.9641879 8.6289465E-02 O.OOOOOOOE+OO 0.2623782 

-0.2163476 -0.4403836 1.527500 -6.2253475E-02 

-6.5925181E-02 0.9641879 0.1993573 -0.1412431 O.OOOOOOOE+OO 

-0.1638800 -0.2929416 1.527500 -0.3346617 

-7.2402343E-02 0.9706651 0.1993573 -0.1408039 -1.1130215E-02 

0.0000000E+00 -0.1716609 1.545834 -0.3346617 

-7.2402343E-02 0.9770384 0.1929839 O.OOOOOOOE+OO -1.1084885E-02 

1.0034929E-03 -0.1716609 1.569711 -0.1939274 

-7.2518051E-02 0.9770384 0.1930996 -2.0943689E-05 O.OOOOOOOE+OO 

1.0032743E-03 -0.1756709 1.569711 -0.1903024 

-7.2518289E-02 0.9770387 0.1930996 -2.0943688E-05 1.0010066E-08 

O.OOOOOOOE+OO -0.1764212 1.569627 -0.1903024 

-7.2518289E-02 0.9770387 0.1930996 O.OOOOOOOE+OO 1.0010066E-08 

-1.3371427E-13 -0.1764212 1.569629 -0.1902815 

-7.2518289E-02 0.9770387 0.1930996 -2.5195758E-21 O.OOOOOOOE+OO 

-1.3371427E-13 -0.1764212 1.569629 -0.1902815 

-7.2518289E-02 0.9770387 0.1930996 -2.5195758E-21 -1.6049783E-34 

O.OOOOOOOE+OO -0.1764212 1.569629 -0.1902815 

O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO O.OOOOOOOE+OO 
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**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 


THE NEW COEFF. OF X SQUARE TERM IS : -7.2518289E-02 
THE NEW COEFF. OF Y SQUARE TERM IS : 0.9770387 
THE NEW COEFF. OF Z SQUARE TERM IS : 0.1930996 
THE NEW COEFF. OF X TERM IS : -0.1764212 
THE NEW COEFF. OF Y TERM IS : 1.569629 

THE NEW COEFF. OF Z TERM IS : -0.1902815 
THE NEW CONSTANT OF PROP. IS : 0.6643000 

**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 

ABCFGHPQR 

0.04714 0.60888 0.44160 -0.86181 0.34007 0.00000 -0.44038 1.13892 -1.01980 

0.04714 0.96419 0.08629 0.00000 0.26238 -0.21635 -0.44038 1.52750 -0.06225 

-0.06593 0.96419 0.19936 -0.14124 0.00000 -0.16388 -0.29294 1.52750 -0.33466 

-0.07240 0.97067 0.19936 -0.14080 -0.01113 0.00000 -0.17166 1.54583 -0.33466 

-0.07240 0.97704 0.19298 0.00000 -0.01108 0.00100 -0.17166 1.56971 -0.19393 

-0.07252 0.97704 0.19310 -0.00002 0.00000 0.00100 -0.17567 1.56971 -0.19030 

-0.07252 0.97704 0.19310 -0.00002 0.00000 0.00000 -0.17642 1.56963 -0.19030 

-0.07252 0.97704 0.19310 0.00000 0.00000 0.00000 -0.17642 1.56963 -0.19028 

-0.07252 0.97704 0.19310 0.00000 0.00000 0.00000 -0.17642 1.56963 -0.19028 

-0.07252 0.97704 0.19310 0.00000 0.00000 0.00000 -0.17642 1.56963 -0.19028 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 

**************************************************** ************ 

**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 
**************************************************** ************ 

Alpha Beta Gamma 


-7.688372 39.50777 40.75703 

-4.519698 5.172774 -1.195895 

2.7384600E-02 7.6535594E-04 1.0796233E-06 
-3.6497606E-12 O.OOOOOOOE+OO O.OOOOOOOE-fOO 

**************************************************** ************ 

**************************************************** ************ 
ALPTOT BETTOT GAMTOT 
-12.18069 44.68131 39.56113 

**************************************************** ************ 
**************************************************** ************ 
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**************************************************** ************ 
**************************************************** ************ 


THE ROTATION MATRIX IS : 
0.8062407 - 0.3102598 - 0.5037009 

- 0.1032203 0.7646025 - 0.6361828 

0.5825129 0.5649087 0.5844286 
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